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Abstract 

If  a  linear  program  (LP)  possesses  a  large  generalized  network  (GN)  subma- 
trix,  this  structure  can  be  exploited  to  decrease  solution  time.  The  problems  of 
finding  maximum  sets  of  GN  constraints  and  finding  maximum  embedded  GN 
submatrices  are  shown  to  be  NP-complete  indicating  thai  reliable,  efficient  solu- 
tion of  these  problems  is  difficult.  Therefore,  efficient  heuristic  algorithms  arc 
developed  for  identifying  such  structure  and  are  tested  on  a  selection  of 
twenty-three  real-world  problems.  The  best  of  four  algorithms  for  identifying  GN 
constraint  sets  finds  a  set  which  is  maximum  in  twelve  cases  and  averages  99.1% 

C"> 

of  maximum.  On  average,  the  GN  constraints  identified  comprise  more  than 
62.3%  of  the  total  constraints  in  these  problems.  The  algorithm  for  identifying 
embedded  GN  submatrices  finds  submatrices  whose  sizes,  rows  plus  columns, 
average  96.8%  of  an  LP  upper  bound.  Over  91.3%  of  the  total  constraint  matrix 
was  identified  as  a  GN  submatrix  in  these  problems,  on  average. 

"The  act  of  being  wise  is  the  act  of  knowing  what  to  overlook." 

William  James   (ca.  1890) 

1.   Introduction 

Large-scale  linear  programming  (LP)  models  frequently  have  sparse 
coefficient  matrices  with  special  structure.  If  special  structure  can  be 
identified,  it  can  often  be  exploited  to  reduce  the  cost  of  solving  the  LP.  "Direct 
factorization,"  e.g.  [13],  maintains  a  partitioning  of  the  rows  and/or  columns  of 
all  simplex  bases.  Computations  are  reduced  with  respect  to  standard  methods 
if  special  structure  can  be  isolated  within  the  partitions.  "Decomposition,"  e.g. 
[14],  splits  a  problem  into  a  master  problem  and  one  or  more  subproblems.  This 
technique  is  most  efficient  when  subproblems  consist  entirely  of  special  struc- 
ture allowing  their  rapid  solution.  The  details  of  these  exploitation  schemes  will 
not  be  discussed  here. 
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Useful  structures  found  embedded  in  a  subset  of  the  rows  and/or  columns 
of  an  LP  constraint  matrix  im  hide  simple  upper  bounds  (at  most  one  nonzero 
element  in  each  row),  generalized  upper  bounds  (GUI!)  (at  most  one  nonzero 
coefficient  in  each  column),  and  networks  (at  most  two  nonzero  elements  in 
each  column).  Varieties  of  embedded  networks  include  the  general  ease,  gen- 
eralized networks  (GN),  generalized  transshipment,  networks  (GT)  (at  most  one 
coefficient  not  equal  to  +1),  and  pure  networks  (NET)  (at  most  one  +1  and  one  -1 
in  each  column). 

Simple  upper  bounds,  GUD  and  NET  structures  have  been  exploited  in  vari- 
ous commercial  and  experimental  optimization  systems,  and  efficient  automatic 
identification  schemes  have  been  developed  to  find  these  structures,  e.g., 
[4,7,8]. 

Recent  research  has  produced  very  efficient  specialized  simplex  algorithms 
for  solving  network  problems  (For  example,  sec  [[}]  for  NET,  |(5]  for  GN,  and 
[6, 11]  for  GT.).  This  research  has,  in  turn,  been  exploited  to  develop  factorized 
optimization  systems  which  solve  general  LP  problems  with  a  set  of  rows  exhibit- 
ing NET  structure  [12],  GN  structure  [IB],  and  GT  structure  [  12,19].  Even  more 
recently,  optimization  systems  have  been  tested  which  use  direct,  factorization 
[19]  or  primal  and/or  dual  decomposition  [14]  to  cxploi!  embedded  GN  struc- 
ture. 

Now  that  software  is  available  to  solve  GN  (and  GT)  problems  [6],  it  is  very 
Likely  that  several  research  groups  will  exploit  GN  in  various  wa\rs  in  the  near 
future.  To  support  this  research,  we  are  interested  in  efficiently  and  automati- 
cally identifying  GN  structure  of  the  following  varieties  in  g<  neral  LP  coefficient 
matrices: 

GNc  A  subset  of  LP  columns  which  are  GN,  or 

GNR  A  subset  of  LP  rows  which  are  GN,  or 

G\:n.e         An  embedded  G\  within  a  subset  of  the  rows  and  columns  of  IP. 

iusc  the  efficiency  of  solving  a  general  LP  with  GN-cxploitin^  methods  is 
enhanced   if  the  GN  structure  is   large,   maximum  GN  structures  are  cur  goal. 
■  the  maximization  problems  described  below. 

Let  li=-\o-jl  be  the  mxn  coefficient  matrix  of  LP,  and  let  H=[/^,-J  be  the 
associated  0-1  incidence  matrix  for  A.  The  three  maximization  problems,  formu- 
lated as  integer  programs,  are 

M(G\:):  max  V  c; 


c      j 


st.        V  h^jCj  <  2    for  all  ;' 

i 

where  c3  is  a  binary  decision  variable  indicating  inclusion  of  column  j  in  GN'c; 
A/(GNR):  maxVr, 

sc         E^'ri-2    for  all; 

i 

7\  e   10,15, 
where  rt  is  a  binary  decision  variable  indicating  inclusion  of  row  i  in  GNH;  and 
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A/(GNR.C):  rnax£rt    +  £  c> 

s.t.  Yj  hijri  +  rrijCj   <  2  +rrij     for  all  j 

i 

rit  Cj  e  JO. lj. 


where  rt  and  c^  are  binary  decision  variables  indicating  respective  inclusion  of 
row  i  and  column  j  in  GNr  q,  and  where  rrij  =  X!  ^-y  —  2-    Note  that  our  definitions 

i 
of  maximum  GN   factorizations   are   expressed  simply  as   the   sum  of  the  rows 
and/or  columns  included. 

Much  work  has  been  done  on  the  development  of  algorithms  to  identify  spe- 
cial substructures  in  LPs.  Previous  work  in  identifying  GUB  subsets  of  con- 
straints is  well  known  [4,7].  Brown  and  Wright  [8]  have  explored  ways  to  identify 
NET  subsets.  Extraction  of  hidden  NET  structure  with  general  linear  transforma- 
tions has  been  discussed  by  Bixby  and  Cunningham  [2]  and  by  Musalem  [20]. 
Identification  of  GN  row  sets  and  other  structures  has  been  proposed  by  Schrage 
[21]. 

The  problems  of  identifying  maximum  GUB  and  NET  constraint  subsets  are 
NP-complete,  and  consequently,  exact  solutions  cannot  be  guaranteed  to  be 
obtained  quickly.  Since  GUB  and  NET  constraints  are  special  cases  of  GN  con- 
straints, it  is  to  be  expected  that  exact  solutions  of  the  GN  identification  prob- 
lems will  also  be  difficult  to  obtain.  We  show  that  the  GN  identification  problems 
are,  in  fact,  NP-complete,  but  also  give  effective  and  reliable  heuristic  algo- 
rithms for  them. 

In  section  2,  the  complexities  of  the  three  maximization  problems  are 
investigated.  A/ (GNr)  and  A/(GNr.c)  are  shown  to  be  difficult  so,  in  section  3, 
efficient  algorithms  are  developed  for  finding  approximate  solutions  to  these 
problems.  Four  specialized  integer  programming  heuristics  are  described  for 
identifying  maximal  GNr  sets.  Two  of  the  algorithms  are  "addition"  heuristics 
which  begin  with  the  empty  CNr  set  and  successively  add  rows  while  maintaining 
feasibility.  The  other  two  algorithms  are  "deletion"  heuristics  which  begin  with 
an  infeasible  GNr  set  and  successively  delete  rows  until  a  feasible  set  is  found. 
Algorithm  GNRC  for  A/(GNr.c)  takes  as  input  the  GNr  set  found  by  any  one  of  the 
the  GNr  heuristics.  Then,  it  successively  adds  rows  which  introduce  the  least 
amount  of  weighted  infeasibility  and  drops  those  columns  where  an  infeasibility 
results.  In  this  way,  a  sequence  of  GNr.c  sets  is  produced  and  the  maximum  of 
these  taken  to  be  the  heuristic  solution  to  A/(GNr.c)-  After  the  algorithms  are 
presented,  computational  experience  is  given  in  section  4. 

2.   Complexity 

In  this  section  we  investigate  the  complexity  of  Af(GNc),  A/ (GNr),  and 
A/ (GNr.c).  A/(GNc)  is  trivially  solvable  in  polynomial  time  by  choosing  all  columns 
with  at  most  two  nonzero  elements  in  them;  consequently,  its  complexity  will  not 
be  discussed  further.   The  other  two  problems  are  more  interesting. 

Following  standard  practice,  A/(GNR)  and  A/(GNRC)  will  be  studied  with 
respect  to  their  associated  decision  problems: 

Z?(GNr):        Does  there  exist  a  set  of  rows  R  in  II  such  that,  for  positive  integer 
k<m, 

|  R  |  >  k     and      £  hy  <  2    for  all  j  ? 
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Z?(GNrc):  Does  there  exist  a  set  of  rows  R  and  columns  C  in  H  such  that,  for 
positive  integer  A;  <m+n, 

\R\  +  \C\&k     and      £/ii;<2    foralljeC? 

ieff 

Of  course,  a  polynomial  algorithm  for  one  of  the  above  decision  problems  would 
imply  a  polynomial  algorithm  for  the  associated  maximization  problem  using, 
say,  a  binary  search  on  the  values  of  k . 

We  consider  the  complexity  of  Z?(GNr  c)  first.  Yannakakis  [24]  investigated 
the  problem  of  finding  the  least  number  of  nodes  which  can  be  deleted  from  a 
bipartite  graph  such  that  the  resulting  induced  subgraph  has  a  particular  pro- 
perty. Restated  in  terms  of  the  decision  problem,  he  gives  the  following 
theorem  on  0-1  matrices  as  a  corollary  of  his  results  on  graphs. 

Theorem  1:  Let  Q  be  any  class  of  0-1  matrices  which  is  closed  under  permuta- 
tion and  deletion  of  rows  and  columns.  Let  H  be  an  raxn  0-1  matrix,  and  let  A; 
be  some  positive  integer,  A:  <  m+n.  Then,  finding  an  m0xnG  submatrix  Hq  of  II 
such  that  Hq  e  Q  and  m0-i-n0  >  A;  is  polynomial  if  the  matrices  of  Q  have  bounded 
rank  and  is  NP-complete  otherwise. 

It  is  assumed  above  that  membership  in  Q  can  be  determined  in  polynomial  time 
for  a  matrix  of  bounded  size  (otherwise,  NP-hardness  would  be  implied). 

This  theorem  is  impressive  in  that  it  handles  the  NP-completeness  question 
for  0-1  matrices  in  a  wholesale  fashion.  The  NP-completeness  of  Z?(GNrc)  follows 
as  a  simple  corollary. 

Corollary  1:   Z?(GNriC)  is  NP-complete. 

Proof:  Let  Q  be  the  class  of  0-1  matrices  with  at  most  two  Is  in  each  column.  Q 
is  obviously  closed  under  permutation  and  deletion  of  rows  and  columns; 
matrices  of  arbitrarily  large  rank  can  be  found  in  Q;  and  membership  in  Q  can 
be  determined  in  polynomial  time.  Z?(GNRC)  for  the  incidence  matrix  II  is 
equivalent  to  searching  for  an  inQxnQ  submatrix  1I0  of  II  such  that  IJjtQ  and 
7n0+n0>A:.   Therefore,  by  Theorem  1,  Z}(GNrC)  is  NP-complete.  ° 

A  0-1  matrix  H  is  represented  as  a  bipartite  graph  with  nodes  on  one  side  of 
the  bipartition  corresponding  to  rows,  nodes  on  the  other  side  of  the  bipartition 
corresponding  to  columns,  and  an  edge  (i,j)  for  each  /iy  =1.  /}(GNr) 
corresponds  to  a  node-deletion  problem  with  deletions  restricted  to  one  side  of 
the  bipartition;  Yannakakis's  results  do  not  directly  apply  since  they  pertain  to 
node  deletions  on  either  side  of  the  bipartition1.  Therefore,  we  use  a  problem- 
specific  proof  to  show  that  Z?(GNr)  is  NP-complete. 

Lemma  1.    £>(GNr)  is  NP-complete. 

IVoof:  For  ease  of  representation,  ZJ(GNr)  will  be  equivalently  stated  in  matrix 

not  ,il 

Z?(GNR):    Does  there  exist  a  binary  m -vector  x  such  that  lx>  A;  and  Hrx<  2? 

ZJ(GNr)  is  obviously  in  NP.  We  show  that  it  is  NP-complete  by  a  transformation 
from  the  "Exact  Cover  by  3-Sets"  problem  [15],  as  specialized  by  Garey  and 
Johnson  [10]. 

D(X3C):  Does  there  exist  a  brnary  p-vector  y  such  that  I3'  =  q  and  Ny  =  1 
where  N  is  a  3gxp,  0-1  matrix  with  exactly  three  Is  in  each  column 
and  at  most  three  Is  in  each  row? 


Bartholdi  [l]  has  addressed  this  topic,  but  his  results  are  incomplete.    For  instance, 
without  additional  restrictions,  his  Theorem  2  would  imply  that  £>(GNC)  is  NP-complete. 
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For  each  row  i  in  N  with  only  one  1  or  two  Is,  augment  N  with  one  or  two  unit 
vector  columns  eit  respectively.  Since  none  of  these  columns  could  be  included 
in  an  exact  cover  of  size  q,  D(X3C)  is  equivalent  to 

D(X3C):       Does  there  exist  a  binary  vector  y'  of  length  p  +Z  such  that  ly*  =  q 

and  (E.N)y'  =  1  where  E  corresponds  to  I  augmenting  columns? 
By  construction  of  D(X3C),  no  set  of  columns  of  cardinality  less  than  q   could 
ever  cover  all  the  rows  exactly  once  let  alone  more  than  once.    Thus,  D(X3C)  is 
equivalent  to  a  "minimum  cover  problem" 

D(UC):         Does  there  exist  a  binary  vector  y"  such  that  iy  <  q  and  (E.N)y'  >  1? 

Let  x=l-y".  Since  each  row  contains  exactly  three  Is,  D(UC)  is  equivalent  to  a 
"maximum  uncover  problem" 

Z?(MUC):       Does    there    exist    a    binary    vector    x   such    that    lx>  p  +  l-q     and 
(E,N)x<2? 

Since  all  above  transformations  are  of  polynomial  complexity,  and  since  D(MUC) 
is  an  instance  of  £>(GN;k),  /J(GNr)  is  NP-complete.  ■ 

3.   Algorithms 

The  complexity  results  of  the  preceding  section  indicate  that  solving 
jW(GNr)  and  A/(GNr.c)  exactly  could  be  very  time-consuming.  Therefore,  heuris- 
tic algorithms  have  been  developed  for  obtaining  approximate  solutions.  We 
describe  the  algorithms  for  M(GNr)  first. 

i/(GNft)  is  an  integer  programming  problem  of  the  form  max  cxs.t.  Axr?b, 
x  binary,  where  all  data  is  nonnegative.  Thus,  integer  programming  heuristics 
seem  appropriate  for  attacking  this  problem.  Two  basic  heuristic  techniques 
exist  for  solving  such  integer  programs  which  we  label  "addition"  heuristics  and 
"deletion"  heuristics.  An  addition  heuristic  begins  with  the  feasible  solution  x=0 
and  successively  sets  to  1  that  variable  Xj  which  which  myopieally  maximizes 
effective  profit.  The  effective  profit  associated  with  x}  is  Cj/tpj,  where  </>j  is  a 
penalty  whose  definition  varies  between  heuristics,  but  which  in  some  way 
reflects  the  units  of  feasibility  used  up  by  setting  x}  to  1.  The  addition  heuristic 
stops  when  when  no  additional  variables  can  be  set  to  1  without  violating  feasibil- 
ity. A  deletion  heuristic  begins  with  the  usually  infeasible  solution  x=  1  and  suc- 
cessively sets  to  0  that  variable  Xj  which  myopieally  minimizes  loss  of  effective 
profit  Cj/(pj.  Here,  pj  is  a  penalty  which  reflects  the  amount  of  infeasibihty 
currently  being  contributed  by  Xj  =  l.  The  deletion  heuristic  stops  when  when  a 
feasible  solution  is  obtained. 

We  have  specialized  two  addition  heuristics  and  two  deletion  heuristics  to 
A/(GNr).  The  addition  heuristics  begin  with  an  empty  GNr  set  and  successively 
add  rows  to  the  set  until  a  maximal  set  is  obtained.  The  deletion  heuristics 
begin  with  an  infeasible  GNr  set  consisting  of  all  the  rows,  and  rows  are  succes- 
sively deleted  until  a  feasible  set  is  obtained.  Since  a  GNr  set  obtained  by  dele- 
tion may  not  be  maximal,  a  second  phase,  an  addition  phase,  is  appended  to 
insure  that  the  set  is  maximal.  To  further  expand  the  GNr  set  found,  it  is  possi- 
ble to  devise  post-maximal  techniques  similar  to  the  2-opt,  3-opt  and  general  k- 
opt  procedures  used  in  traveling  salesman  heuristics,  e.g.,  [16,17].  Application 
of  such  techniques  was  unwarranted,  however,  since  computational  results  in 
section  4  show  that  excellent  approximate  solutions  were  obtained  using  the 
basic  addition  and  deletion  heuristics. 

The  addition  heuristics  are  described  by  Algorithm  GNRa,  with  variations 
"Greedy"  and  "Toyoda."  The  effective  profit  associated  with  adding  row  ?'.  to  the 
GNr  set  is  1/ RPi  where  RPi  is  a  row  penalty  derived  from  the  current  nonmaxi- 
mal   solution,  the   nonzero   elements   in  the  row  and  feasibility   requirements. 
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Thus,  at  each  step  of  the  algorithm,  the  row  with  the  smallest  penalty  is  added 
to  the  GNr  set.  Feasibility  is  maintained  by  setting  to  infinity  the  row  penalty  of 
any  row  whose  addition  would  cause  an  infeasibility.  In  the  Greedy  variation, 
RPi  equals  the  number  of  nonzero  elements  in  the  row  if  the  penalty  is  finite. 
The  Toyoda  variation  is  a  modification  of  an  integer  programming  heuristic 
developed  by  Toyoda  [23].  In  this  heuristic,  the  finite  row  penalty  RPi  is  based 
not  only  on  the  number  of  nonzero  elements  in  the  row,  but  also  on  how  close  to 
feasibility  limits  addition  of  the  row  would  bring  the  current  solution. 

The  deletion  heuristics  are  described  by  Algorithm  GNRd,  with  variations 
"Dobson"  and  "Scnju  &  Toyoda."  In  this  algorithm,  each  row  has  a  penalty  RPl 
which,  roughly  speaking,  indicates  how  much  infeasibility  the  row  is  contribut- 
ing. 1/  RPi  is  the  loss  in  effective  profit  if  row  i  is  removed  from  the  GN#  set. 
Thus,  this  algorithm  successively  deletes  rows  with  maximum  penalty  to  minim- 
ize the  loss  of  effective  profit. 

Dobson  [9]  analyzes  and  gives  worst-case  performance  guarantees  for  an 
addition  heuristic  for  integer  programs  of  the  form  min  ex,  s.t.  Ax^b,  0<x<u,  x 
integer,  where  all  data  is  nonnegative.  By  simple  substitution  of  variables,  how- 
ever, the  Dobson  heuristic  may  be  interpreted  as  a  deletion  heuristic  for  prob- 
lems in  the  form  of  //(GNr).  At  each  deletion  step  of  this  heuristic,  RPt  is  the 
number  of  nonzero  elements  in  row  i  which  are  contributing  to  an  infeasibility. 
If  ra0  is  the  optimal  solution  to  4/ (GNr)  and  mD  is  the  heuristic  solution  obtained 
by      deletion       only,       Dobson's      worst-case      bound      on      performance      is 

(m  —  mj))/ (ttl  —  m0)  <  ^   \/k    where  d   is   the  maximum  number  of  nonzero 

fc  =  i 
elements  in  any  row.  This  is  the  only  performance  guarantee  known  for  any  of 
the  heuristics  implemented  in  this  paper.  Unfortunately,  the  upper  bound  on 
ttlq  this  yields  is  rather  weak  in  practice  (See  Table  3.).  Any  addition  heuristic 
may  be  used  as  a  second  phase  for  a  deletion  heuristic,  but  for  the  Dobson  dele- 
tion heuristic,  we  chose  the  greedy  addition  heuristic  as  the  second  phase  since 
the  definition  of  RPt  is  consistent  between  the  two  phases. 

The  second  variant  of  GNRd  is  a  specialization  of  the  heuristic  devised  by 
Senju  and  Toyoda  [22]  which  those  authors  label  an  "effective  gradient  method." 
For  A/(GNR),  l\T  maps  the  set  of  feasible  r  values  into  the  n. -dimensional  hyper- 
cube  whose  sides  are  of  length  2.  At  every  step  of  the  algorithm,  given  current 
infeasible  solution  r,  RPi=(liTr-2)*-bl,  where  the  jfth  element  of  (Hrr-2)+  is 

n 

max  $0,  £]Aiyri— 2j.    RP\  may  be  interpreted  as  the  length  of  the  projection  of 

*=i 
the  vector  h*  onto  the  shortest  vector  extending  from  the  point  \\Tr  outside  of 

the  hypcrcube  to  the  boundary  of  the  hypercube.  The  modified  Toyoda  addition 
heuristic  is  used  as  the  second  phase  of  this  heuristic. 

The  two  algorithms  GNRa  and  GNRd,  with  their  variations,  are  outlined  as 
follows: 


-7- 


Algorithm  GNRa 

Input:  The  LP  coefficient  matrix  A. 

Output:  A  set  of  row  indices  7R  corresponding  to  the  largest  GNr  set  found  in 
A 

Comment:  The  basic  algorithm  is  the  "Greedy"  addition  heuristic.  The  modified 
"Toyoda"  heuristic  is  obtained  by  substituting  the  statement  in 
square  brackets  for  its  predecessor. 

Step  0.      "Initialization" 
Initialize: 

(a)  I  =  <f>  and  /'  =  $1,2 raj. 

(b)  For  each  column  j,  a  column  bound 

CBj  =  2 

Comment:    CBj  is  the  number  of  elements  column  j  may  con- 
tain. 

(c)  For  each  ie/',  a  row  penalty 

RPt  =    £  1 

Step  1.      "Row  Addition" 

Let  RP  -  RPS  be  the  smallest  row  penalty  (corresponding  to  rowse/'). 
If  Rp  <  oo  then 

(a)  Move  s  from  /'  to  /. 

(b)  For  each  column  j  such  that  aSJ-  /  0, 

(i)     Let  CBj  -  CB}  -  1. 

(ii)   If    CBj  =0    then    for    each    i?s    such    that    ay  *  0,    let 
RPt  -  °°. 

[{ii)   For    each    i^s    such    that    a^-  ^  0,    if    CBj-l    then    let 
/?/*  =  RPi  +  1,  else  let  RP,  =  «. 

(c)  Repeat  Step  1. 

Step  2.      "Termination" 

Print  /R  =  /  and  STOP. 

End  of  Algorithm  GNRa 
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Algorithm  GNRd 

Input:  The  LP  coefficient  matrix  A. 

Output:  A  set  of  row  indices  /r  corresponding  to  the  largest  GNr  set  found  in 
A 

Comment:  The  basic  algorithm  is  the  "Dobson"  heuristic.  The  "Senju  and  Toy- 
oda"  heuristic  is  obtained  by  substituting  the  statements  in  square 
brackets  for  their  predecessors. 

Step  0.      "Initialization" 
Initialize: 

(a)  /  =  \\,2 mj  and/'  =  0. 

(b)  For  each  column  j,  a  column  penalty 

CPj  =  (    S    1  )  "  2. 

ie/ 

Comment:    CPj  is  the  number  of  "excess"  elements  in  column 


(c)     For  each  ie/,  a  row  penalty 

RPi  =    S    l. 

CPj>0 

Comment:  RP,  is  number  of  units  of  infeasibility  which  row  i 
is  currently  contributing. 

(c)    For  each  iel,  a  row  penalty 

RPi  =    S    cpi- 

By  1*0 

CPj>0 

Comment:    RPi  is  the  sum  of  excess  elements  in  columns  with 
_       a  nonzero  entry  in   row  i. 

Step  1.      "Row  Deletion" 

Let  RP  =  RPi  be  the  largest  row  penalty  (corresponding  to  row  £e/). 
If  PP  >  0  then 

(a)  Move  Z  from  /  to  /'. 

(b)  For  each  column  j  such  that  atj  ^  0 

(i)     If    CPj  =  1    [If    CPj  >  0]    then,    for    each   i*s    such   that 
Oy ,  /  0,  let  /?Fi  =  RPi  -  1. 

(ii)   Let  CPj  =  CPj  -  1. 

(c)  Repeat  Step  1. 

Step  2.      "Row  Addition  Penalties" 

For  each  ie/',  compute  a  row  penalty 


RP, 


RP,  =, 


£  1     if  CP;  <  0  for  all  Og  /  0 

°°  otherwise. 

2    (CPj  +  3)     if  CPj  <  0  for  all  o^  *  0 

Oy*0 

00  otherwise. 
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Step  3.      "Row  Addition" 

Let  RP  -  RPS  be  the  smallest  row  penalty  (corresponding  to  row  seP). 
URP<™,  then 

(a)  Move  s  from  /'  to  /. 

(b)  For  each  j  such  that  asj  *  0,  let  CPj  -  CPj  +  1. 

(c)  Go  to  Step  2. 

Step  4.      "Termination" 

Print  /r  =  /  and  STOP. 

End  of  Algorithm  GNRd 

The  execution  times  of  the  above  algorithms  and  the  other  algorithms 
described  in  this  paper  are  quite  short  if  proper  data  structures  are  used.  The 
initial  computation  of  the  row  and  column  penalties  can  be  made  very  quickly  if 
the  nonzero  entries  in  each  row  and  column  are  stored  in  a  linked  list.  Column 
penalties  are  then  updated  in  a  single  pass  of  a  row.  Because  of  sparsity,  row 
penalties  can  usually  be  updated  in  passes  through  just  a  few  columns. 
Efficiency  is  further  improved  if  row  and  column  partitions  are  maintained  with 
an  indirect  address  array  which  allows  contiguous  access.  Associated  with  this 
mapping  array,  a  second  array  expresses  the  inverse  map  to  speed  updating. 

An  easily  computable  upper  bound  on  A/(GNr),  denoted  UBr,  is  useful  for 
checking  the  efficacy  of  the  above  algorithms.  Algorithm  UBR  is  designed  for 
this  purpose.  Let  A]  and  A2  be  a  partition  of  the  rows  of  A  and  let  z ,  z,  and  zz  be 
the  solutions  to  M(GNr)  on  A  At  and  A2,  respectively.  If  UBX  is  any  valid  upper 
bound  on  M(GNr)  for  A],  then 

z<z,  +  z2<  UBi  +  z2. 

Algorithm  UBR  iteratively  applies  the  above  statement,  computing  the  simple 
bound  UB\  and  letting  A  =  Ag  after  each  iteration.  This  is  repeated  until  all 
columns  of  Aq  have  at  most  two  nonzero  elements  in  them  at  which  point  z2  is 
equal  to  the  number  of  rows  in  A?.  UBr  is  then  given  by  the  sum  of  the  UBX 
upper  bounds  found  at  each  iteration  plus  zz  found  at  the  last  iteration.  At  each 
iteration,  A  is  partitioned  with  respect  to  that  column  j  having  the  maximum 
number  of  nonzero  entries.  A2  is  all  rows  of  A  with  Oy  #  0  and  VBX  =  zx  =  2  since 
column  jf  has  only  nonzero  elements  in  Aj. 


Algorithm  UBR 

Input:        The  LP  coefficient  matrix  A. 
Output:     A  value  UB^,  an  upper  bound  on  |  /r| 

Step  0.      "Initialization" 
Initialize: 


(a)  /  =  1 1,2 to  J,  and  UBR  -  0. 

(b)  For  each  column  ,7,  a  column  count 

cq  =  (  2  1 ). 


id 


Step  1.      "Iterative  Partitioning" 

Let  CC  -  CQ  be  the  largest  column  count  (corresponding  to  column  I). 
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If  CC  >  2  then 

(a)  Let  UBR  =  UBR  +-  2. 

(b)  For  each  iel  such  that  a^  #  0, 

(i)     Delete  i  from  /. 

(ii)   For  each  j  such  that  Oy  ^  0,  update  column  count  letting 
CCj  -  CCj  -  1. 

(c)  Repeat  Step  1. 

Step  2.      "Termination" 

Print  UBR  =  UBR  +  |  7 1  and  STOP. 

End  of  Algorithm  UBR 

Algorithm  GNRC,  the  heuristic  for  7/(GNr  c)>  is  outlined  next.  Any  one  of  the 
integer  programming  heuristics  described  for  A/(GNr)  could  be  applied  to  this 
problem.  However,  these  algorithms  will  normally  give  only  a  single  answer  to 
the  problem;  our  algorithm  allows  the  exploration  of  a  complete  trajectory  of 
maximal  GNrc  sets  beginning  with  GNr  and  ending  with  GN^.  Our  algorithm 
begins  with  the  set  of  rows  7R  found  in  Algorithm  GNRa  or  GNRd  and  repeatedly 
attempts  to  expand  this  set  by  deleting  columns,  always  saving  the  largest  GNrC 
set  found.  This  approach  was  suggested  by  manual  analysis  of  several  problems 
for  which  the  GNr  set  is  limited  by  a  few  key  complicating  columns.  Deleting 
these  columns  produced  a  much  larger  embedded  GNrC  set,  and  motivated 
development  of  a  new  factorization  LP  code  which  effectively  exploits  GNRC 
structure  [19]. 


Algorithm  GNRC 

Input         The  LP  coefficient  matrix  A  and  a  GNr  set  /r,    |/r|   <  m,  e.g.,  Ir  from 
Algorithms  GNRa  or  GNRd. 

Output:     A  set  of  row  indices  7r_c  and  a  set  of  column  indices  Jrc  corresponding 
to  the  largest  GNr  c  structure  found  in  A. 

Siep  0.      "Initialization" 
Initialize: 

(a)  /  =  FR,    V  =  J  1.2 in\  -I,    J  =  J  1,2 n\,   7RC  =  /,   and 

/RC  =  J  ■ 

Comment:    7  and  J  are  the  current  sets  of  row  and  column 
indices  while  7r.c  and  JRq  store  the  best  sets  found. 

(b)  For  each  column  jtJ,  a  column  penalty 

cp> ■>  =  (  E  i )  -  2- 

Comment:    These  column  penalties  remain  as  an  artifact  of 
Algorithm  GNRd  and  can  be  defined  as  input. 

(c)  For  each  iel',  a  row  cost 

RCi  =    2    1. 

Comment:     RCt    is   the   number   of   columns   which   must  be 
deleted  if  row  i  is  added  to  I . 
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Step  1.      "Column  Deletion" 

Let  RC  =  RCS  be  the  smallest  row  cost  (corresponding  to  rowse/'). 

(a)  For  eachje/  such  that  cu^  *  0, 

(i)     Let  CPj  =  CPj  +  1. 

(ii)   If  CPj  =  1  then  delete  j   from  /  and  for  each  ieP  such 
that  dij  it  0,  update  row  costs  letting  RQ  =  RCt  -  1. 

(b)  Move  s  from  /'  to  I. 

Step  2.      "Row-inclusion  Penalties" 

For  each  ie/',  compute  a  row  penalty 

£    (CPj  +  1)     if  CPj  <  0  for  all  By  ^  0 


RPi  = 


°°  otherwise 


Step  3.      "Row  Addition" 

Let  RP  -  RPS  be  the  smallest  row  penalty  (corresponding  to  row  se/'). 
If  i?P<0  then 

(a)  Move  s  from  /'  to  /. 

(b)  For  each  jzj  such  that  as,-  ^  0 

(i)     Let  CPj  =  CPj  +  1. 

(ii)   If  CP;  =  0  then  for  each  ie/'  such  that  Oy  ^  0,   update 
row  costs  letting  RQ  =  RCX  +•  1. 

(c)  Go  to  Step  2. 

Step  4.      "Incumbent  Test" 

If  |  /  |  +  |  J"  |  >  |  Tr.c  I  +  |  J\i,c\  then  let  7RC  =  /  and  J"R  c  =  J. 
Step  5.      "Termination" 

If  |  /  |  <  77i ,  then  go  to  Step  1.    Otherwise,  print  /rc,  Jr.c  and  STOP. 
End  of  Algorithm  GNRC 

A  stronger  test,  allowing  preemptive  termination,  is  possible  at  Step  5:  If 
1 1 1  <  m  and  m.  +-  |/|  >  |/R_cl  +  I^R.cl-  However,  the  weaker  test  permits  the 
exploration  of  a  complete  trajectory  for  GNR  c  as  discussed  above. 

Along  the  lines  of  W?R,  an  easily  computed  upper  bound  on  A/(GNRc). 
denoted  UB^q,  was  developed  to  check  the  accuracy  of  GNRC.  Partition  A  as  fol- 
lows: 


A  = 


An  I  A12 


A,i  I  A 


-d 


Let  z,  2,,  and  222  be  the  solutions  to  yW(GNR)  on  A  An  and  A22.  respectively,  and 
let  UBn  be  any  simple  upper  bound  on  A/(GNRC)  for  An.   Then 

2  <  zn  +  zZ2<  UBn  +  222. 

Algorithm  UBRC  computes  UBRC  by  iteratively  applying  the  above  statement, 
computing  the  simple  bound  UBn  and  letting  A=  Aq2  after  each  iteration.  This 
is  repeated  until  all  columns  of  A22  have  at  most  two  nonzero  elements  in  them 
at  which  point  22  is  equal  to  the  number  of  rows  plus  the  number  of  columns  in 
Ag2.  UBjiC  is  then  given  by  the  sum  of  the  UBX  upper  bounds  found  at  each 
iteration  plus  222  found  at  the  last  iteration.  If  An  is  selected  such  that  it  con- 
sists of  single  column  and  three  rows,  all  with  nonzero  elements,  then 
UBu-zn-3.    Computational  experience  has  indicated  an  effective  rule   for 
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selecting  the  partition:  Among  all  columns  in  A  having  at  least  3  nonzeros,  select 
that  column  having  the  minimum  number  of  nonzeros,  and  within  that  column 
select  the  first  three  rows  with  nonzeros  in  them.  If  k  partitions  are  carried  out 
before  /I22  becomes  a  GN  matrix,  it  follows  that: 

UBR,C  =  3*  +  |/ga|  +  |/a8|  =  3k+(\I\-3k)+(\J\-k)  =  |/| +  |J"| -A:. 

The  last  equality  is  used  in  computing  UBRc- 

Algorithm  UBRC 

Input:        The  LP  coefficient  matrix  A 

Output:     A  value  UBRC,  an  upper  bound  on  |  7RC|  +  I  «/r,cI  • 

Step  0.      "Initialization" 
Initialize: 

(a)  I  =  [1,2,  .  .  .  ,mj,  and  UBR_C  =  |/|  +  |  J|. 

(b)  For  each  column  j,  a  column  count 

ccj  =  (  2  1 ). 

id 

Step  1.      "Iterative  Partitioning" 

Let     CC  =  CCS     be     the     smallest     column     count     greater     than     2 

(corresponding  to  column  s). 

If  no  such  column  exists,  go  to  step  2.    Else, 

(a)  Let  UBR,C  =  UBR.C  -  1. 

(b)  For  exactly  three  ie/  such  that  a^  /  0, 

(i)     Delete  i  from  I. 

(ii)    For  each  j  such  that  a^  ^  0,  update  column  count  letting 
CC)  =  cq  -  1. 

(c)  Repeat  Step  1. 

Step  2.      "Termination" 

Print  UBRC  and  STOP. 

End  of  Algorithm  UBR 

4.    Computational  Experience 

The  algorithms  described  in  section  3  have  been  implemented  in  FORTRAN, 
using  the  X-System  [5j  as  the  host  optimization  package.  Table  1  identifies 
twenty-three  LP  and  mixed  integer  programming  (MIP)  problems  which  have 
been  collected  from  various  sources  over  the  years.  Some  of  these  models  are 
very  well  known,  e.g.,  Dantzig's  PILOT  and  the  U.S.  Department  of  Energy's  PAD 
and  PIES,  and  most  of  them  were  sent  to  us  because  of  their  difficulty,  solution 
expense,  or  outright  solution  failure  on  commercial  optimization  systems.  Table 
1  shows  problem  dimensions  excluding  right-hand  sides  and  objective  functions. 
Computation  times  displayed  in  Tables  2-4  are  compute-seconds,  accurate  to 
the  precision  shown,  for  FORTRAN  IV  H  (Extended)  with  0ptimize(2),  run  on  IBM 
3033AP  under  VM/CMS. 
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Table  1.   LP /MIP  Problem  Set 

Problem 

Constraints 

Variables 

Nonzero  Els. 

Model 

AIR 

170 

3,040 

6,023 

Physical  Distribution 

ALUMINUM 

4,045 

6,805 

27,917 

Econometric  Production  &  Distribution 

COAL 

170 

3,753 

7,506 

National  Energy  Planning 

CUBIC  1 

G57 

3,074 

15,894 

Combinatorics  Problem 

CUBIC2 

2,689 

11,905 

63,361 

Bigger  Combinatorics  Problem 

CUPS 

360 

618 

1,341 

Production  Scheduling 

ELEC 

784 

2,800 

8,462 

Energy  Production  &  Consumption 

FERT 

605 

9,024 

40,484 

Production  &  Distribution 

FOAM 

999 

4,020 

13,083 

Production  Scheduling 

FOOD 

4,010 

14,409 

23,332 

Production,  Distribution  &  Inventory  Planning 

GAS 

788 

5,541 

31,020 

Production  Scheduling 

JCAP 

2,486 

3, 849 

9,510 

Production  &  Shipment  Scheduling 

LANG 

1,235 

1,425 

22,026 

Equipment  &  Manpower  Scheduling 

NETTING 

80 

190 

388 

International  Currency  Exchange 

ODSAS 

4,647 

4,995 

30,832 

Manpower  Planning 

PAD 

694 

3,297 

15,541 

Energy  Allocation,  Distribution  &  Consumption 

PAPER 

2,868 

5,318 

23,746 

Econometric  National  Production 

PIES 

662 

3,011 

13,376 

Energy  Production  &  Consumption 

PILOT 

974 

2,172 

12,927 

Energy  Development  Planning 

REFINE 

5,220 

5,994 

40,207 

Oil  Refinery  Model 

STEEL 

831 

1,276 

9,608 

Econometric  Production  &  Distribution 

TRUCK 

220 

4,752 

30,074 

Fleet  Dispatch  (Set  Cover) 

WADDING 

2,991 

15,001 

82,706 

Muiticommod'ty  Prod.  &  Distribution  Planning 

Algorithms  GNRa  and  GNRd  were  used  to  identify  GNr  rows  with  Algorithm 
UBR  used  to  give  an  upper  bound  on  the  total  number  of  such  rows.  To  check 
accuracy,  we  attempted,  within  budget  limitations,  to  solve  exactly  the  integer 
linear  programs  for  A/(GNr)  in  those  cases  where  \Ih\<UBr  ■  (V/e  were  success- 
ful in  all  but  one  case,  as  seen.  Times  for  solving  the  ILPs  averaged  214.1 
seconds  for  those  problems  solved.)  Results  for  GNRa  and  GNRd,  given  in  Table  2, 
are  (a)  the  size  of  the  optimal  GNr  set  found  by  the  ILP,  (b)  the  size  of  this  set  as 
a  percentage  of  total  problem  rows  m,  (c)  the  size  of  the  GNr  set  found  by  GNR, 
(d)  the  size  of  this  set  as  a  percentage  of  the  ILP  optimum,  and  (e)  the  time 
required  by  the  algorithm.  For  GNRd,  the  column  labeled  |  /r|  uses  the  notation 
a:b  where  a  is  |/r|  and  b  is  the  number  of  rows  in  /r  which  were  gained  in  the 
addition  phase  of  the  heuristic.  Problems  are  weighted  equally  in  computing 
average  percentages  in  the  "Totals"  row  of  the  table.  Times  listed  do  not  include 
input  or  output. 

All  GNR  variants  perform  quite  well.  The  addition  phase  in  GNRd  did  not 
often  contribute  a  significant  fraction  of  the  GN  rows  found,  but  the  additional 
rows  found  helped  make  both  GNRd  variants  slightly  better  than  either  of  the 
GNRa  variants.  The  best  algorithm  on  this  problem  set,  GNRd  (Senju  &  Toyoda), 
finds  an  average  of  99.1%  of  the  maximum  GNr  set  on  those  problems  which  we 
can  solve  exactly.  The  GNr  sets  average  62.3%  of  the  total  problem  rows  on 
these  same  problems.  GNR  computation  times  are  nominal  compared  with 
actual  solution  times  of  the  seminal  LPs  and  MIPs. 
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Results  for  UBR,  given  in  Table  3,  include  (a)  the  size  of  the  optimal  GNr  set, 
(b)  the  upper  bound,  (c)  the  upper  bound  as  a  percentage  of  the  ILP  optimum, 
and  (d)  the  time  required  to  find  the  upper  bound.  For  comparison,  we  include 
(e)  Dobson's  upper  bound  labeled  "  UBDr"  and  (f)  that  bound  as  a  percentage  of 
the  ILP  optimum.  Table  3  also  displays  some  properties  of  GNr  as  found  by 
GNRd,  Senju  &  Toyoda.  These  properties  include  (g),  the  number  of  disjoint 
embedded  GN  components,  (h)  the  largest  and  smallest  components,  (i)  the 
number  of  null  columns,  and  (j)  the  number  of  singleton  columns.  These  proper- 
ties are  of  interest  since  the  structure  of  the  embedded  generalized  network 
affects  the;  solution  techniques  used  in  an  LP  factorization.  For  example,  com- 
ponents consisting  of  single  rows  may  be  handled  most  efficiently  without  utiliz- 
ing a  complete  generalized  network  code. 


Tabic  3    CNR  Features 

#(GNR) 

AJg 

onLhri 

UBR 

Dobson 

Bound 

Embedded  G 

Mp  Components 

Problem 

ILP  Opt 

<JBR 

%  Opt 

Time 

UBDR 

%  Cpt. 

Total 

Largest 
(m+n) 

Smallest 

Null 
CcIj 

Sins 

AIR 

170 

170 

1C0 

.0 

170 

1C0 

1 

17043.040 

0 

57 

ALUMINUM 

2.193 

2.214 

1C0? 

1.7 

3.798 

172  8 

1  15 

1.1  I  Hh  3.431 

1  i  I 

0 

1 .03  1 

COAL 

170 

170 

100 

0 
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1C0 
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1?C<  3.753 

C 

0 

CUBIC  1 
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324 

103  0 
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190.7 
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10CH.719 

H-2 
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G12 

CUBIC2 
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1C5  1 
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10G  1 
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G.034 
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1C0G 

.0 

091 
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1    5+7.1  17 

1  +  4 

532 
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H  4 
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PAD 
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.0 

4  0?  0 

3 

00  i  1.001 

8  +  33 
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PAPER 
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■1 

2730 
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225H.GC1 

1  1  1 

C75 

1 .70  1 
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.0 

571 

193  3 

35 

11C*  -1,915 
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9.0G 
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101  G 
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■  ..I.  ■ 
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17?  C 

95 
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1C5 
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.1 
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1  +  18 

1.7C3 
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2,211 

o  noo 

ICC  5 

.2 

2800 

123  6 

3 

909  i -1.1  CO 

1+1 
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5. CO  2 

UBr  is  surprisingly  tight,  averaging  101.4%  cf  the  true  maximum,  and  com- 
putation times  are  nominal.  Dobson's  bound  is  poor,  tng  203.1%  of  th.2 
true  maximum.  The  GN  components  found  usually  consist  of  a  few  large  com- 
ponents and  numerous  small  components. 
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Table  4  gives  the  results  obtained  by  Algorithm  GNRC  and  Algorithm  UBRC. 
Since  no  ILP  optimum  is  known  for  ^/(GNrc)  in  most  cases,  the  items  displayed 
differ  from  those  items  displayed  in  Tables  2  and  3.  The  results  reported  for 
Algorithm  GNRC  are  (a)  the  size  of  the  of  the  GNrc  structure  found,  (b)  the  time 
in  seconds  required  to  find  the  structure  excluding  input  and  output,  (c)  the  size 
of  the  GNrc  as  a  percentage  of  the  total  constraint  matrix,  and  (d)  the  percen- 
tage of  of  the  LP  upper  bound  (UBLPrq)  achieved  by  the  algorithm.  The  results 
reported  for  Algorithm  UBRC  are  (e)  |/r,cI  +  I^r.cI  as  a  percentage  of  UBr£<  and 
(f)  the  time  required  to  obtain  UBr,c-  For  comparison,  the  last  two  columns  of 
the  table  give  the  total  number  of  rows  and  columns  obtained  for  the  GNr  and 
GNc  problems.  These  are  the  sizes  of  the  embedded  GN  submatrices  when  res- 
tricted to  row  submatrices  and  column  submatrices,  respectively.  Each  prob- 
lem is  weighted  equally  to  compute  average  percentages  in  the  "Totals"  row. 


Table  4. 

GN'R.C  Results 

Algorithm  GNRC 

Alg.  UBRC 

GNR 

GNC 

Problem 

KrcI  +  I«/r.cI 

Time 

%  (m+n) 

%  UBLPRC 

%  UBRC 

Time 

\IR\+n 

\JC\+rn 

AIR 

3,210 

.0 

100 

100 

100 

.0 

3,210 

3,210 

ALUMINUM 

9,02? 

13.6 

83.2 

91.7 

91.7 

2.3 

8,980 

5,508 

COAL 

3,923 

.0 

100 

100 

100 

.0 

3,923 

3,923 

CUBIC  1 

3,365 

.6 

90.2 

99.4 

94.8 

.4 

3,365 

659 

CUBIC2 

13,096 

11.1 

89.7 

99.5 

94.7 

6.4 

13,096 

2,690 

CUPS 

951 

.0 

97.2 

100 

99.7 

.0 

951 

713 

ELEC 

3,322 

.3 

92.7 

99.0 

98.3 

.2 

3,320 

1,042 

FERT 

9,596 

.3 

99. 7 

100 

99.9 

.2 

9,595 

2,362 

FOAM 

4,971 

.1 

99.0 

100 

99.7 

.1 

4,871 

1,044 

FOOD 

18,137 

8 

98.5 

99.5 

99.4 

.1 

18,125 

17,860 

GAS 

5,920 

5.4 

93.5 

94.9 

94.5 

.2 

5,614 

848 

JCAP 

5,822 

5.5 

91.9 

97.7 

99.8 

.2 

4,651 

5,718 

LANG 

2,139 

1.1 

60.4 

97.8 

90.2 

.2 

2,133 

1,905 

[NG 

262 

.0 

93.9 

97.6 

100 

.2 

262 

256 

ODSAS 

7,556 

40.0 

78.4 

78.0 

86.1 

1.2 

6,470 

5,094 

PAD 

3,621 

3.9 

00.7 

98.8 

95.3 

.3 

3,419 

2,416 

PAPER 

7,388 

4.6 

89.9 

95.9 

96.2 

.9 

7,179 

4,905 

3,313 

.9 

90.2 

99.5 

94.8 

.2 

3,299 

2,241 

OT 

2,645 

1.4 

84.1 

95.7 

91.6 

.2 

2,634 

1,567 

[NE 

9,326 

19.3 

63.2 

93.8 

92.4 

2.3 

9,104 

7,729 

1.700 

.9 

60.7 

91.5 

69.7 

.2 

1,605 

1,131 

4,822 

97.0 

NA 

90.3 

,o 

4,872 

270 

)ING 

17,209 

6.3 

95.6 

99.7 

97.8 

1.0 

17,209 

14,451 

Totals 

141,321 

118.1 

91.3% 

OS.  8% 

95.6 

17.1 

133,232 

87,582 

Note:    NA  indicates  not  available. 


GNRC  performs  very  well,  also.  The  algorithm  finds  a  GNrc  structure  whose 
size  overages  91.3%  of  the  size  of  the  total  constraint  matrix.  The  size  of  the 
structure  averages  96.8%  of  the  LP  upper  bound  on  those  problems  for  which  the 
bound  was  obtained.  (Times  to  obtain  the  LP  bound  averaged  315.6  seconds.) 
With  respect  to  UBRC,  the  GNRc  set  found  averages  95.6%.  Thus,  the  upper 
bound  provided  by  algorithm  UBRC  is  only  slightly  weaker,  on  average,  than  the 
LP  upper  bound.   In  addition,  UBRC  has  more  than  a  400  to  1  computational 
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speed  advantage  over  the  LP  upper  bound  making  it  very  attractive. 

Additional  computational  studies  have  been  performed  to  investigate  the 
structures  which  GNR  and  GNRC  obtain.  Figure  1  summarizes  this  work  for 
ELEC.  JCAP,  PAD,  PIES  and  PILOT.  The  outer  rectangle  represents,  to  scale,  the 
constraint  matrix  for  each  problem.  The  area  above  the  dashed  line  represents 
the  GNr  set  found  by  GNRd,  Senju  &  Toyoda.  Within  this  area  are  indicated  the 
connected  components  found  by  a  simple  connectivity  algorithm.  As  indicated 
previously  in  Table  3,  a  few  large  components  are  typically  found  together  with 
numerous  small  components.  The  area  to  the  left  of  the  vertical  line  represents 
the  GNc  set.  The  irregular  lines  trace  the  trajectories  of  the  GNrc  structures 
found  by  GNRC,  ranging  from  GNR  on  the  right  to  GNc  at  the  lower  left.  From  any 
point  on  this  trajectory,  all  rows  and  columns  above  and  to  the  left  form  a  GN 
set.   The  circle  indicates  the  largest  GNrc  structure  found  on  this  trajectory. 

5.  Conclusion 

Although  GNC  identification  is  easy,  GNR  and  GNR  c  identification  is  theoreti- 
cally difficult.  However,  maximal,  and  often  optimal  GNr  and  GNRC  substruc- 
tures can  be  found  in  an  LP  constraint  matrix  using  the  heuristic  algorithms 
developed  here.  In  some  problems,  large  GNr  structures  can  be  found,  while  in 
other  problems,  it  is  necessary  to  remove  some  columns  to  find  a  large  embed- 
ded GNRC  structure.  Since  execution  time  is  modest  for  heuristic  GN 
identification,  our  algorithms  can  be  applied  as  a  matter  of  course  in  general 
LPs  to  seek  GN  substructures.  Evidence  from  the  problem  set  indicates  that 
this  is  well-advised  if  a  GN-exploiting  method  is  available:  No  members  of  the 
problem  set  were  known,  a  priori,  to  contain  significant  GN  structure  and  yet,  in 
several  cases,  GN  structure  was  predominant. 
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Figure  1.    Embedded  Generalized  Networks 
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Figure  1.    (continued) 
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